Lag analysis at a state-level

Use most correlated shift time to build linear models

\[y = \beta_{0} + \beta_{1}(\sum I(A_i = 1)) + \beta_{2}x_1 + \beta_{3}x_2\]

# select best lag number based on intial exploration
case.shifted.days.spearman <- 37 # based spearman correlation
doc.visit.shifted.days.spearman <- 62 # based spearman correlation for doctor visit
cum.death.shifted.days.spearman <- 55
death.shifted.days.spearman <- 55

case.vec <- c("confirmed_7dav_incidence_num", "confirmed_7dav_cumulative", "confirmed_7dav_cumulative_prop")
doc.vec <- c("smoothed_cli", "smoothed_adj_cli")
cum.death.vec <- c("deaths_7dav_cumulative_num")
death.vec <- c("deaths_7dav_incidence_num")

# Make two copies for shifting the data
selected.df <- intervention_mobility_case
factored_data <- intervention_mobility_case
# Change the data by shifting the covariates
# Shift the case count column vector by the specified shift time
factored_data <- shiftDays(selected.df, factored_data, case.shifted.days.spearman, case.vec)

# Shift the  doctor column vector by the specified shift time
factored_data <-shiftDays(selected.df, factored_data, doc.visit.shifted.days.spearman, doc.vec)
  
# Shift the death count column vector by the specified shift time
factored_data <-shiftDays(selected.df, factored_data, cum.death.shifted.days.spearman, cum.death.vec)

# Shift the cum death count column vector by the specified shift time
factored_data <-shiftDays(selected.df, factored_data, death.shifted.days.spearman, death.vec)
  
# We specifically look at emergency declaration
factored_data%>% 
  mutate(EmergDec.duration = cumsum(EmergDec)) %$% 
  lm(full_time_work_prop ~ EmergDec.duration + smoothed_cli+smoothed_adj_cli ) %>% 
  summary()
## 
## Call:
## lm(formula = full_time_work_prop ~ EmergDec.duration + smoothed_cli + 
##     smoothed_adj_cli)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0176227 -0.0056012  0.0008011  0.0053044  0.0172350 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        3.429e-02  1.584e-03  21.648  < 2e-16 ***
## EmergDec.duration  8.692e-05  2.250e-05   3.864 0.000153 ***
## smoothed_cli      -4.964e-03  8.586e-04  -5.782 2.95e-08 ***
## smoothed_adj_cli   5.073e-03  1.094e-03   4.635 6.58e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.007426 on 192 degrees of freedom
##   (458 observations deleted due to missingness)
## Multiple R-squared:   0.29,  Adjusted R-squared:  0.2789 
## F-statistic: 26.13 on 3 and 192 DF,  p-value: 3.202e-14
# We specifically look at emergency declaration
factored_data%>% 
  mutate(EmergDec.duration = cumsum(EmergDec)) %$% 
  lm(full_time_work_prop ~ EmergDec.duration + smoothed_cli+smoothed_adj_cli+confirmed_7dav_cumulative_prop+deaths_7dav_cumulative_num) %>% 
  summary()
## 
## Call:
## lm(formula = full_time_work_prop ~ EmergDec.duration + smoothed_cli + 
##     smoothed_adj_cli + confirmed_7dav_cumulative_prop + deaths_7dav_cumulative_num)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0187374 -0.0062279  0.0006184  0.0050205  0.0176432 
## 
## Coefficients:
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     4.089e-02  4.814e-03   8.494 7.72e-15 ***
## EmergDec.duration              -1.150e-05  7.275e-05  -0.158    0.875    
## smoothed_cli                   -5.080e-03  8.776e-04  -5.789 3.14e-08 ***
## smoothed_adj_cli                5.691e-03  1.211e-03   4.701 5.18e-06 ***
## confirmed_7dav_cumulative_prop -4.713e-06  6.172e-06  -0.764    0.446    
## deaths_7dav_cumulative_num      2.012e-06  1.672e-06   1.203    0.230    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.007515 on 178 degrees of freedom
##   (470 observations deleted due to missingness)
## Multiple R-squared:  0.2717, Adjusted R-squared:  0.2512 
## F-statistic: 13.28 on 5 and 178 DF,  p-value: 5.333e-11

Model with shifting the covariates (weekends included)

Model without shifting the covariates (weekends included)

intervention.lm <- intervention_mobility_case %>% 
  mutate(EmergDec.duration = cumsum(EmergDec)) 

lm.fit.no.lag <- lm(full_time_work_prop ~ EmergDec.duration + smoothed_cli+smoothed_adj_cli, data =intervention.lm
) 

intervention.lm$predlm <- c(rep(NA, nrow(intervention.lm) - length(predict(lm.fit.no.lag))), predict(lm.fit.no.lag))

intervention.lm%>% 
  mutate(policy.duration = cumsum(EmergDec), EmergDeclaration = as.factor(EmergDec)) %>% 
  ggplot(aes(x = time_value, y = full_time_work_prop, color = EmergDeclaration)) +
  geom_point() + 
  geom_line(aes(x = time_value, y = predlm, colour="fitted value"), size = 1)+
  labs(title = "Covariates selected WITHOUT most correlated number of shift")

Weekend effects

We suspect that the mobility signal is lower than usual during the weekend.

intervention_mobility_case$weekday <- weekdays(as.Date(intervention_mobility_case$time_value)) 

p <- ggplot(intervention_mobility_case, aes(x=weekday, y=full_time_work_prop)) + 
  geom_boxplot()
p

Re-plot the regression line after dropping weekends

Covariates shifted

## 
## Call:
## lm(formula = full_time_work_prop ~ EmergDec.duration + smoothed_cli + 
##     smoothed_adj_cli)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0202756 -0.0028614 -0.0000421  0.0029460  0.0125748 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        3.867e-02  1.190e-03  32.491  < 2e-16 ***
## EmergDec.duration  7.145e-05  1.688e-05   4.232 4.24e-05 ***
## smoothed_cli       8.542e-04  7.575e-04   1.128    0.261    
## smoothed_adj_cli  -8.345e-04  8.877e-04  -0.940    0.349    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.0047 on 136 degrees of freedom
##   (328 observations deleted due to missingness)
## Multiple R-squared:  0.4607, Adjusted R-squared:  0.4488 
## F-statistic: 38.72 on 3 and 136 DF,  p-value: < 2.2e-16

Covariates not shifted

## 
## Call:
## lm(formula = full_time_work_prop ~ policy.duration + smoothed_cli + 
##     smoothed_adj_cli)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0307231 -0.0062396 -0.0004375  0.0054900  0.0241455 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       5.588e-02  1.290e-03  43.338  < 2e-16 ***
## policy.duration   4.248e-05  1.466e-05   2.898 0.004225 ** 
## smoothed_cli      2.076e-03  1.036e-03   2.003 0.046712 *  
## smoothed_adj_cli -5.100e-03  1.324e-03  -3.852 0.000163 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.009018 on 180 degrees of freedom
##   (284 observations deleted due to missingness)
## Multiple R-squared:  0.1653, Adjusted R-squared:  0.1514 
## F-statistic: 11.88 on 3 and 180 DF,  p-value: 3.896e-07